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We have performed ab initio tensile tests of bulk Al along different tensile axes, as well as perpen- 
dicular to different grain boundaries to determine mechanical properties such as interface energy, 
work of separation and theoretical strength. We show that all the different investigated geometries 
exhibit energy-displacement curves that can be brought into coincidence in the spirit of the well 
known universal binding energy relationship curve. This simplifies significantly the calculation of 
ab initio tensile strengths for the whole parameter space of grain boundaries. 



I. INTRODUCTION 



The prediction of the mechanical properties of poly- 
crystals requires knowledge about the mechanical prop- 
erties of all interfaces, i.e. grain boundaries, in their 
microstructure. So far, most mesoscale models of 
microstructure-property relationships as used in contin- 
uum simulations of deformation and fracture make rather 
simple assumptions for the variation in grain boundary 
properties with the boundary geometryH21. With ab ini- 
tio electronic structure calculations it is possible to de- 
termine the cohesive energy, elastic modulus, sliding bar- 
rier, and theoretical strength of interfaces accurately and 
quantitatively. The considerable computational effort, in 
comparison e.g. to atomistic simulations employing em- 
pirical potentials, is easily feasible with modern comput- 
ers if the investigations are restricted to grain boundary 
structures based on coincidence-site lattices. The ab ini- 
tio approach is desirable to avoid problems of transfer- 
ability of empirical potentials, usually fitted to reproduce 
equilibrium properties, to non-equilibrium processes such 
as failure. It is indispensable whenever a phenomenon is 
controlled by the electronic structure. This is the case 
in systems with directional bonds, or when studying the 
influence of chemical composition, and alloying effects. 
Nevertheless, sampling the five parameter space given 
by the geometric degrees of freedom of the grain bound- 
ary (rotation axis and angle, and the grain boundary 
normaP}) by ab initio calculations remains a challenge. 
Models exist that relate the energies of certain subsets 
of this space to the grain boundary geometry. Most 
well known is the dislocation model for small angle tilt 
grain boundaries based on the picture of Read and Shock- 
ley, which can be extended empirically to large-angle tilt 
grain boundaries- . If the energy of a grain boundary 
can be related to its geometric parameters via such mod- 
els, the description of the energy hypersurface is signif- 
icantly simplified. However, so far it seems that there 
is no such correlation which would be valid in the com- 
plete parameter space^l The complexity of the problem 
is increased by the fact that for use in mesoscale mod- 
els we are not only looking for a function that describes 
the grain boundary energy as function of misorientation, 



but also for its first and second derivatives with respect 
to a displacement from the equilibrium volume, i.e. the 
maximum stress and the elastic modulus. Nevertheless, 
the situation is not hopeless, as the analytic function of 
this energy-displacement curve is the same for any grain 
boundary geometry, if obtained under the same loading 
conditions. The demonstration of this fact, a universal 
binding behaviour in the spirit of Rose's universal bind- 
ing energy relationship (UBER/^, is subject matter of 
this paper. 

In section|Tl]we describe our computational procedure, 
including an explanation of grain boundary nomenclature 
and details on different ways to perform il ab initio" ten- 
sile tests. In section Hill we review the UBER and its im- 
plications. In the results' section, sec. |IV| we give details 
about the grain boundary structures after a full optimiza- 
tion of the microscopic degrees of freedom ( IV A ) and the 



corresponding energies ( IV B ) . The results of the tensile 
tests are presented in section |IV C| The universal elastic 
behaviour under tensile load of all systems investigated 



is demonstrated in section IV D and the relationship be- 
tween energies and strength is discussed in section |IV E| 
We summarize our insights in section [V] 



II. TECHNICAL DETAILS 

For a given orientation of the tensile axis, we con- 
structed supercells for bulk, surface, and grain boundary 
calculations of the same size and shape. In other words, 
starting from a bulk supercell containing N atomic lay- 
ers, half of the planes were replaced by vacuum to create 
a surface slab, or half of the planes were replaced by the 
same number, but with a misorientation, to create a grain 
boundary structure. 

In detail, supercells for tensile tests were constructed 
for Al bulk, such that the z-axis is oriented along the 
[111], [112], [113], and [114] direction, i.e. defining (111), 
(112), (113) and (114) as the cleavage plane. In addi- 
tion we constructed special grain boundaries containing 
these planes as grain boundary planes. These were the 
S3 (111) [111] 60° twist grain boundary, the £3 (112) 
[110] 109° symmetrical tilt grain boundary (STGB), the 
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Sll (113) [110] 129° STGB, and the £9 (114) [110] 141° 
STGB. In this nomenclature the use of £ indicates that 
for the chosen misorientation a periodic superstructure 
can be found, the so-called coincidence site lattice (CSL). 
Therefore these grain boundaries are also called special 
grain boundaries. The value of £ is the volume of a unit 
of this CSL divided by the volume of the cubic Al unit 
cell, i.e. it is a measure for the periodicity of the grain 
boundary. The grain boundary plane is given in round 
brackets, the direction of the axis of misorientation in 
rectangular ones. All tilt grain boundaries considered 
here are symmetrical, which means that the grain bound- 
ary plane divides the misorientation angle in two equal 
parts. In other words, it represents a mirror plane. Note 
that the (111) [111] 60° twist grain boundary can also be 
expressed as (111) [110] 70.5° tilt grain boundary, but we 
prefer to refer to it as twist grain boundary to emphasize 
its close-packed atomic structure. 

To calculate the total energy of the above mentioned 
supercells we performed density-functional-theory (DFT) 
calculations of total energy and electronic structure em- 
ploying the ABINIT open source code^. The exchange 
correlation effects were treated in the local density ap- 
proximation and electron-ion interactions were modeled 
via a Trouller-Martins type norm-conserving pseudopo- 
tential for Al. Convergency with respect to fc-point den- 
sity, plane wave energy cut-off, and system size was tested 
for the (111) surface energy. The surface energy is given 
by 
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where Ef^. is the total energy of the surface slab and 
-Efot the energy of the bulk supercell containing twice 
the number of atomic layers. A is the surface area. The 
plane wave cut-off was varied between 12 and 20 Ha in 
steps of 2 Ha. To test convergency with respect to cell 
size, we chose bulk (surface) supercells containing 18 (9), 
24 (12), or 30 (15) (111) planes. The fc-point meshes 
employed were of the Monkhorst-Pack type, using 2x2x1, 
4x4x1, 8x8x2, and 12x12x2 fc-points. With a plane wave 
cut-off of 16 Ha, the 8x8x2 Monkhorst-Pack mesh, and 
a minimum surface slab thickness of nine atomic layers, 
the surface energy was converged within an accuracy of 
MO- 4 Ha/atom (w 2.5-10~ 3 eV/atom). 

For orientations different from z parallel to [111] it was 
ensured that the distance between the interfaces was at 
least as large as in the [111] cells. The supercell shapes 
were made commensurate with that of the [111] cell, 
which also enables the use of a commensurate fc-point 
mesh. 

To obtain accurate interface energies all microscopic 
degrees of freedom were optimized and the positions of 
the atoms were relaxed until the remaining forces were 
smaller than 5-10" 5 Ha/a.u. («3-10^ 3 eV/A ). After- 
wards, the excess interplanar spacing of the grain bound- 
aries can be calculated as half the difference (due to the 
periodic boundary conditions) between the relaxed (Do) 



and initial (Dhki) supercell length perpendicular to the 
interface, 
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The grain boundary energy is given by 
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where E^ is the total energy of the grain boundary 
supercell and E^ k the energy of the bulk supercell con- 
taining the same number of crystallographic layers in the 
corresponding orientation. A is the interface area. From 
the total energy of the surface slabs, also the work of 
separation of a grain boundary was calculated according 
to 
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Note that the work of separation for the perfect bulk cor- 
responds to twice the surface energy as defined in eq. ([I]) . 

After a fit of the energy-displacement curves, the ten- 
sile strength can be calculated as the slope in the inflec- 
tion point: 
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where A is the displacement from the equilibrium in- 
terplanar distance. Initially, the tensile tests were per- 
formed in three different ways: (a) by a simple scaling 
of the supercell dimensions along the tensile axis, while 
leaving the internal coordinates fixed. Atomic relaxations 
then lead to a homogeneous strain distribution in bulk 
supercells, and a characteristic distribution of strain in 
any cell containing a defect (e.g. a grain boundary), (b) 
By performing rigid grain shifts (rgs). Here the spacing 
between two blocks of atoms is increased only at a de- 
fined cleavage plane (due to the use of periodic boundary 
conditions at two defined cleavage planes per supercell). 
Within the blocks, the interplanar distance corresponds 
to the equilibrium bulk value. This way we can model 
ideally brittle cleavage under loading mode I. (c) By do- 
ing the latter and relaxing the atomic positions at each 
shift, while keeping the total elongation of the supercell 
fixed. Again, this corresponds to a mode I cleavage pro- 
cess, but now elastic energy is released due to atomic 
relaxations. 



III. UNIVERSAL BEHAVIOR AND SCALING 
LENGTHS 

Rose et alP postulated and demonstrated that the 
binding energies Et of metals have a universal form of 
the kind 



E h (d) = \E° h \g(a) 



(6) 
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where d is the interatomic distance, or, in the case of 
interface energies the interplanar spacing. \E^\ is the 
binding energy at equilibrium volume (— — W sep ) and a 
is the rescaled displacement 



The characteristic length scale I depends on the curvature 
of the energy- volume curve at the minimum, i.e., on the 
elastic modulus as follows, 

If rescaled in this manner, all energy-volume curves coin- 
cide, i.e. they have the same functional form 17(a). This 
phenomenon was observed for adhesiorP and cohesion^ 
of metals, as well as chemisorption on metal surfaces^. 
More recently, Hayes et al. have shown that even 
the energy-displacement curves of non metallic systems 
(AI2O3 and Si) can be brought into coincidence with 
those of metals (AYp^l. This universal behavior of chem- 
ically different systems means that we can determine the 
cohesive behavior (i.e. theoretical strength and critical 
displacement) of any material from three parameters, E^, 
do, and E h (do), once the functional form g(a) is known. 
To describe hydrostatic volume expansion / compression 
of simple metals, Rose determined g(a) to be 

g(a) = -(l + a + 0.05a 3 )e~ a . (9) 

Hayes et al. used an asymptotic approximation, i.e. a 
simple quadratic function introduced by Nguyen et alP^l 
that scales with the system size to fit results of uniax- 
ial computational tensile tests. Originally applied to the 
results of homogeneous strain, the approach was gener- 
alized by Hayes et al. to fit the results of tensile tests 
in the form of rgs + atomic relaxation, thus taking into 
account surface relaxations^. The simple function and 
the scalability of the results make this approach attrac- 
tive for the calculation of traction-separation laws used in 
continuum models. However, if the theoretical strength 
shall be calculated independently of the system size, a 
function displaying an inflection point is to be preferred. 
To fit the results of our tensile tests in the form of rigid 
grain shifts, we used 

g(a) = -(l + a)e- a , (10) 

a function also used by Rose, to represent the results of 
displacing metal-metal interfaces of different metals. 

IV. RESULTS 

A. Grain boundary structures 

After construction according to the macroscopic pa- 
rameters the microscopic degrees of freedom of the grain 





do [A] 


h [A] 


lr[k] 


(111) bulk 




0.561 


1.652 


(112) bulk 




0.653 


1.599 


(113) bulk 




0.655 


1.695 


(114) bulk 


- 


0.665 


1.657 


E3(lll)[lll]60° twist GB 


0.014 


0.615 


1.745 


E3(ll2)[110]109° STGB 


0.389 


0.621 


1.736 


Ell(ll3)[110]129° STGB 


0.311 


0.620 


1.508 


E9(ll4)[110]141° STGB 


0.185 


0.605 


1.514 



TABLE I: Excess volume at the different grain boundaries, 
represented by do, and scaling lengths for brittle and relaxed 
cleavage, and l r , for bulk and grain boundary planes (ex- 
planation in the text). 



boundaries were optimized by performing rigid grain 
shifts perpendicular as well as parallel to the interface, 
followed by relaxation of the atomic positions. All grain 
boundaries exhibit excess volume at the interface, cal- 
culated according to equation (j2|). These expansions da 
are localized in the vicinity of the grain boundary, but 
not necessarily confined only to the first crystallographic 
plane at the interface. The expansions are summarized in 
table |l] The smallest expansion is observed for the (111) 
twist grain boundary, which is the most dense cleavage 
plane after the (111) bulk plane. The stable translation 
state of the Ell STGB parallel to the interface is the ini- 
tial, mirror-symmetric one. In the case of the E3 STGB 
a shift of 0.5 times the interplanar spacing along the tilt 
axis [110] breaks this mirror-symmetry and produces a 
structure about 150 mJ/m 2 lower in energy than the 
mirror-symmetric one. The E9 STGB initially contained 
two atomic columns at unphysically small distance at the 
grain boundary. Relaxing the structure (rgs + atomic 
relaxations) led to major re-arrangements of the atomic 
positions at the interface. Removing instead one atomic 
column at the grain boundary and relaxing the struc- 
ture (rgs + atomic relaxation) led to an interface energy 
which is lower by 0.067 Ha/atom (w 1.8 eV/atom). The 
resulting structure at the grain boundary is in excellent 
agreement with experimental observations^. 



B. Interface energies and work of separation 

The resulting interface energies, i.e. surface and grain 
boundary energies, and the work of separation for the 
systems investigated are shown in table |Hj The grain 
boundary energies follow the known trend for the ener- 
gies of symmetric [110] tilt grain boundaries in fee metals 
as function of the misorientation angle, with e nergy c usps 
representing the E3 twist and the Ell STGrP^. Re- 
sults of previous DFT investigations were found for the 
E3 twisiPland the Ell STGB^and agree well with our 
results. The work of separation for the grain boundaries 
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(111) bulk 


0.987 


1.973 


12.6 


(112) bulk 


1.109 


2.218 


12.8 


(113) bulk 


1.093 


2.187 


12.7 


(114) bulk 


1.145 


2.290 


13.2 


E3(lll)[lll]60° twist GB 


0.048 


1.924 


11.5 


E3(112)[110]109° STGB 


0.393 


1.818 


6.8 


Ell(113)[110]129° STGB 


0.171 


2.007 


12.0 


E9(114)[110]141° STGB 


0.486 


1.779 


12.1 



TABLE II: Interface energy, work of separation, and theoret- 
ical strength of Al bulk and grain boundaries. 




A [A] 

FIG. 1: Results of ab initio tensile tests along the [111] di- 
rection in bulk Al. Black circles represent the results of rgs, 
open circles that of rgs followed by relaxation of the atomic 
positions, and squares that of a homogeneous straining of the 
supercell. 



does not follow the same trend, as it also depends on the 
corresponding surface energy. The lowest surface energy 
in Al is that of the close packed (111) plane. Thus, W sep 
for the stable twin is lower that that of the Ell STGB, 
corresponding to the second- most favorable [110] grain 
boundary orientation in foe metals. 



C. Energy-displacement curves 

Starting from the optimized structures, tensile tests 
were performed as explained in section [TTJ Figure [T] shows 
the results for bulk Al stretched and compressed along 
the [111] direction. Rigid grain shifts without any atomic 
relaxation produce a rather steep curve in a shape that 
can be fitted well by equation (10). The energy asymp- 



totically approaches the unrelaxed (111) surface energy, 
which is slightly higher than the relaxed one taken as 
the reference energy level in fig. [T] After relaxing the 
atomic structure at each shift, the atomic positions and 
the corresponding energies coincide with the results ob- 
tained by stretching the supercell homogeneously. A de- 
viation is observed above a relative displacement of the 
two bulk slabs of 2A. In the case of the rgs, the atoms 
relax to a configuration where the strain is localized at 
the defined cleavage plane, whereas in case of the homo- 
geneous strain the crystallographic planes have the same 
distance. A continuation of the homogeneous tensile test 
would lead to N free standing crystallographic planes, if 
N is the number of layers in the supercell. 

If the tensile test is performed with a supercell con- 
taining a grain boundary, the strain is expected to lo- 
calize at this defect also in the case of a homogeneous 
elongation of the supercell. This is shown in figure [2] 
for the £11(113)[110]129° symmetrical tilt grain bound- 
ary. Again, the black circles represent the results of 
the rigid grain shift calculations without atomic relax- 
ations. The open circles mark the energies after atomic 
relaxation. The diamonds indicate the energies obtained 
after a successive scaling of the total length of the su- 
percell, i.e. starting with the equilibrium structure the 
relaxed reduced coordinates of the atoms are taken as 
input coordinates of the next strain state. This ensures 
a continuous deformation path up to the point where 
the energy crosses the reference level, the energy of the 
relaxed (113) surface. At this point the energy drops 
and then starts to increase again. An analysis of the 
atomic structure at this discontinuity shows that the two 
grains have partially debonded and are only connected 
via chains of atoms containing the initial coincidence site. 
Due to this coincidence site a further elongation of the 
cell again leads to a free standing crystallographic plane, 
in the center between two (113) surfaces. Thus, while 
such a drag calculation is interesting to investigate fail- 
ure modes and critical displacements, it is not suitable to 
calculate the work of separation or the tensile strength. 
We can conclude that for the latter a defined point of 
failure is needed in every loading scheme. 



D. Scaling lengths for different geometries 



As mentioned in section |TTl| a universal behavior of ma- 
terials under strain has been shown for bulk metals and 
ceramics, and for coherent metal-metal interfaces. In our 
investigation, we stick to one metal, Al, but probe differ- 
ent orientations of the tensile axis for bulk tensile tests, 
as well as different grain boundary geometries, where the 
tensile axis is perpendicular to the interface. As can be 
seen in figure [3] the results of rgs calculations scale per- 
fectly. The parameters resulting from the fit are given in 
tabic [n] (Eh = —Wsep)- For the case of rgs followed 
by atomic relaxations, the function given in equation 
( 10 1 does not describe the results well. The fit could 
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FIG. 2: Results of ab initio tensile tests at the Ell STGB. 
Black circles represent the results of rgs calculations without 
atomic relaxations, and open circles mark the energies after 
atomic relaxation. Diamonds are the energies obtained via a 
"drag" calculation (see text). The reference level is the energy 
of the relaxed (113) surface. 
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FIG. 3: Rgs energy displacement relationships re-scaled ac- 
cording to eqs. [6] to [8] using eq. |10| For details of the twist 
and STGB geometries see text. 



be improved by using a polynomial including higher or- 
der terms (>3). However, for the rescaling procedure the 
exact function g(a) actually does not have to be known. 
Instead we determined the curvature at minimum en- 
ergy by assuming a quadratic function close to the mini- 
mum, i.e. making a harmonic approximation. The results 
of rgs calculations followed by atomic relaxations were 
rescaled with the curvature thus obtained and with the 
binding energy. The rescaled curves are shown in figure 
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FIG. 4: Energy displacement relationships with relaxed 
atomic positions, re-scaled according to eqs.[5]to[8] using the 
harmonic approximation for g(a). For details of the twist and 
STGB geometries see text. 



[4] Apart from a few outliers the agreement is satisfac- 
tory. Deviations occur mainly at displacements where 
the structure becomes unstable and where we can think 
of a crack forming. Thus, whether investigating relaxed 
or unrelaxed cleavage of different geometries, all we have 
to know about a system to describe its elastic response 
upon tensile load are the binding energy in the specified 
cleavage plane and the elastic constant (resp. E (A)) for 
the tensile direction. Especially for the grain boundaries 
we can conclude that once we have determined a fitting 
function 17(a) we can use it for the complete parameter 
space spanned by the degrees of freedom of the interfaces. 

The scaling lengths calculated to make the different 
energy-displacement curves coincide are given in table 
|lj For ideally brittle cleavage (rgs without atomic relax- 
ation), the results for bulk Al agree very well with the 
empirical results of Rose et al.^, 0.66A. This value was 
obtained by using the experimental results of Simmons 
and WanglSl for the elastic constant of Al along the [111] 
direction and the surface energy of Tyson and Miller^. 
The latter represents an average over different crystal- 
lographic planes, thus we expect our value to be more 
exact. Ab initio calculations for bulk Al have been car- 
ried out by Lazar et alP^. For tensile tests along the [111] 
direction they obtain 4 = 0.54 A and l r = 2.4 A. While 
there is excellent agreement between the scaling lengths 
for brittle cleavage, our result for the relaxed cleavage, 
also shown in table |TJ is about 30% lower. This is due 
to the fact that, although the energy-displacement curve 
is fitted to a quadratic function in both cases, the defi- 
nition of l r is different. Lazar and Podloucky^ni define l r 
as the critical opening, at which the elastic energy E(A) 
equals the work of separation, comparable to the discon- 
tinuity in our drag results, see figure [l] This means their 
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quadratic fit is performed in the spirit of the asymptotic 
approximation of Nguyen and Ortiz^, however it is in- 
dependent of the system size. In our case, where the 
quadratic fit is strictly confined to the region around the 
energy minimum, l r has no such meaning, but is simply 
a measure for the stiffness of the cleavage plane, as in the 
work of Rose et al. 



there is no simple relation between energy and strength, 
and (7th has to be calculated explicitly. 

Fortunately, it is not necessary to calculate full energy- 
displacement curves for the complete parameter space of 
grain boundaries, because, as we were able to show, these 
exhibit a universal behavior. Thus, after having deter- 
mined the analytical function g(a) once, we only need to 



E. Tensile strength 

From the derivatives of the energy-displacement 
curves, the theoretical tensile strength of the bulk phase 
along a given tensile axis and of grain boundaries along 
an axis perpendicular to the grain boundary plane, can be 
calculated according to equation ^ . For direct compar- 
ison of bulk and grain boundary properties, we restrict 
ourselves to the results of rgs calculations without re- 
laxation of the atomic positions. Thus, the strain is fully 
localized between the cleavage planes in both cases. A re- 
laxation would lead to different strain distributions in the 
supercell, as discussed in detail in section [TTJ The results 
are summarized in table [TTJ The theoretical strength of a 
grain boundary is generally lower than that of bulk lattice 
planes in the corresponding orientation. Furthermore, in 
the bulk, the trend in theoretical strength follows the one 
in the work of separation: the higher W scp , the higher 
CTth ■ For the grain boundaries, however, although the 
same overall trend is still visible, the relationship is not 
so straight forward. 



SUMMARY AND CONCLUSIONS 
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FIG. 5: Mechanical properties of the different grain bound- 
aries and bulk supercells as function of the misorientation 
from (110). Energies are represented by circles, and stresses 
by squares. Grain boundary properties are displayed in white, 
and bulk and surface properties in black. For more explana- 
tions see text. 



We have performed ab initio tensile tests of bulk Al 
along different tensile axes, as well as perpendicular to 
different grain boundaries. It has been discussed that 
in order to simulate a physically meaningful de-cohesion 
process, a plane of cleavage has to be defined in all sys- 
tems, also those containing defects as grain boundaries. 
This was done by performing the tensile tests by means 
of rigid grain shifts, followed or not by relaxation of the 
atomic positions at each shift. The mechanical proper- 
ties of the investigated systems are summarized again in 
figure [5} The grain boundary energies agree well with the 
empirical extension of the Read-Shockley picture, show- 
ing cusps at the £3 (111) and the Ell (113) orientations. 
The most stable defect judging from the interface energy 
is the (111) twist grain boundary. However, this is not 
the grain boundary with the highest strength, demon- 
strating that for defect structures as grain boundaries, 



determine the minimum energy structure, corresponding 
surface energies, and the curvature at the minimum. By 
this simplification a use of ab initio tensile strengths in 
continuum models of polycrystals is getting within reach. 
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